
#include <stdlib.h>
#include <stdio.h>
#include <gmp.h>
#include <time.h>


void borwein(unsigned, mpf_t, unsigned int *);

int main(int argc, char *argv[])
{

	unsigned prec_bits = 35000000;
	time_t start, end;
	double dif = 0.0;
	unsigned int i;
	mpf_t pi;
	mpf_set_default_prec(prec_bits);
	mpf_init(pi);

	time(&start);
	borwein(prec_bits,pi,&i);
	time(&end);

	dif = difftime(end,start);

	printf("Time: %2lf\n",dif);
	printf("Iterations: %u\n",i);
    printf("PI = \n");
    mpf_out_str(stdout, 10, 10000001, pi);
    printf("\n");
    mpf_clear(pi);
	return 0;
}


void borwein(unsigned precisionBits, mpf_t pi, unsigned int *i)
{
	mpf_t a;
	mpf_t y;
	mpf_t temp,temp2,temp3;

	mpf_init(a);
	mpf_init(y);
	mpf_init(temp);
	mpf_init(temp2);
	mpf_init(temp3);

	//sqrt(2)
	mpf_sqrt_ui(temp,2);
	
	//a
	mpf_mul_ui(a,temp,4);
	mpf_ui_sub(a,6,a);

	//y
	mpf_sub_ui(y,temp,1);


	for(*i=0;*i < 11;(*i)++)
	{
		mpf_pow_ui(temp,y,4);
		mpf_ui_sub(temp,1,temp);
		mpf_sqrt(temp,temp);
		mpf_sqrt(temp,temp);

		mpf_add_ui(temp2,temp,1);
		mpf_ui_sub(temp,1,temp);

		mpf_div(y,temp,temp2);

		mpf_add_ui(temp,y,1);
		mpf_pow_ui(temp2,temp,4);
		mpf_mul(temp2,a,temp2);

		mpf_pow_ui(temp3,y,2);
		mpf_add(temp,temp,temp3);
		mpf_mul(temp,y,temp);
		mpf_mul_2exp(temp,temp,2*(*i) + 3);
		mpf_sub(a,temp2,temp);
	}

	mpf_ui_div(pi,1,a);
	mpf_clear(a);
	mpf_clear(y);
	mpf_clear(temp);
	mpf_clear(temp2);
	mpf_clear(temp3);
}








